/*
    Spherical Grid Volume Approximation: An algorithm for approximating the 
    volume of a compound object by keeping track of its members' positions in space.
    Copyright (C) 2009, Rodrigo Castano, FuDePAN

    This file is part of Spherical Grid Volume Approximation.

    Spherical Grid Volume Approximation is free software: you can redistribute 
    it and/or modify it under the terms of the GNU General Public License as published 
    by the Free Software Foundation, either version 3 of the License, or
    (at your option) any later version.

    Spherical Grid Volume Approximation is distributed in the hope that it will be useful,
    but WITHOUT ANY WARRANTY; without even the implied warranty of
    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
    GNU General Public License for more details.

    You should have received a copy of the GNU General Public License
    along with this program.  If not, see <http://www.gnu.org/licenses/>.
*/



#ifndef GRILLADO_H
#define GRILLADO_H

#include <cmath>
#include <iostream>
#include <stdexcept>
#include <list>
#include <utility>
#include <algorithm>
#include "tmp_algorithms.h"

using namespace std;

typedef int CustomId;
typedef list< CustomId > CustomIdContainer;
typedef list< CustomId > CustomIdIterator;
//typedef pair < CustomIdContainer , unsigned int > SphereInfo;

typedef float Coord;
typedef float Volume;
typedef float Length;
static const Length PI = atan(1.0) * 4.0f;

class esferaId {
	public:
		esferaId(GridCoord a, GridCoord b, GridCoord c);
		esferaId(const esferaId &s);
		inline GridCoord get_m() const;
		inline GridCoord get_n() const;
		inline GridCoord get_z() const;
	protected:
		GridCoord v;
		GridCoord w;
		GridCoord z;
};
inline esferaId::esferaId(GridCoord input_v, GridCoord input_w, GridCoord input_z) : v(input_v), w(input_w), z(input_z) {}

inline esferaId::esferaId(const esferaId &s) : v(s.get_m()), w(s.get_n()), z(s.get_z()) {}

inline GridCoord esferaId::get_m() const {
	return v;
}
inline GridCoord esferaId::get_n() const {
	return w;
}
inline GridCoord esferaId::get_z() const {
	return z;
}

template <typename T>
class Grillado
{


public:
	//~ struct CustomIdNonExistent : public domain_error {
		//~ CustomIdNonExistent()
			//~ : domain_error( "The CustomId provided does not belong to the indicated sphere." ) {}

	//~ };
	
	template <int X, int Y, int Z, size_t X_MAX, size_t Y_MAX, size_t Z_MAX >
	template <class T_MATRIZ>
	friend unsigned int generic_step<X, Y, Z, X_MAX, Y_MAX, Z_MAX>::apply(T_MATRIZ *** matriz, GridCoord coord_x, GridCoord coord_y, GridCoord coord_z);
	
	struct gridPositionException : public domain_error {
		gridPositionException()
			: domain_error( "The desired position is already empty" ) {}

	};
	struct sizeParamException : public invalid_argument {
		sizeParamException()
			: invalid_argument( "M<3 || N<3 || Z<3" ) {}

	};

	struct radiusOrDistanceParamException : public invalid_argument {
		radiusOrDistanceParamException()
			: invalid_argument( "Either 1) sqrt(2)*radio >= dist: there may be diagonal intersections between spheres (spheres which differ in more than one coordinate) or 2) dist>=2*R: there are no intersections between spheres." ) {}

	};
	
	Grillado(size_t a, size_t b, size_t c, Length R, Length D) throw 
			(radiusOrDistanceParamException, sizeParamException, bad_alloc);

	~Grillado() ;

	
	void info_grillado() const;

	//inline volElementHandler add_element(const volGridAtom &a) ;

	
	
	//nonvolElementHandler add_element( const nonvolGridAtom &a);
	
  //void get_customIDs_container(Coord x, Coord y, Coord z, CustomIdContainer &it) const ;
	
	
	//inline void remove_element(const volElementHandler &id) throw (gridPositionException);
	
	
	//inline void remove_element(const nonvolElementHandler &id) ;
		
	// Devuelve el volumen parcial.
	inline Volume obtener_vol_parcial() const;
	void reset();
	
protected:
	const size_t v;
	const size_t w;
	const size_t z;
	const Length r;
	const Length d;
	const Volume vol_esfera;
	const Volume vol_int;
	
	T ***matriz;
	//SphereInfo ***matriz;	
	unsigned int intersecciones;
	unsigned int esferas;

	inline static Volume calc_vol(Length radio);
	
	inline static Volume calc_int(Length radio, Length dist);
	
	//void add_volumetric_element(const volElementHandler &id);
	
	//inline void add_non_volumetric_element(const nonvolElementHandler &id);
	
	inline unsigned int show_quantity(T &sphere_data) const;
	//inline unsigned int show_quantity(GridCoord x, GridCoord y, GridCoord z) const;
	inline void reset_specific_sphere(T &sphere_data) const; 
	void reducir_vol_parcial(GridCoord coord_x, GridCoord coord_y, GridCoord coord_z) ;
	void aumentar_vol_parcial(GridCoord coord_x, GridCoord coord_y, GridCoord coord_z) ;
	void calcular_coord(Coord coord_x, Coord coord_y, Coord coord_z, 
			   GridCoord &x, GridCoord &y, GridCoord &z) const;
	unsigned int calcular_intersecciones(GridCoord coord_x, GridCoord coord_y, GridCoord coord_z) const;
	//inline volElementHandler get_vol_element_handler(Coord x, Coord y, Coord z, CustomId id) const;
	//inline nonvolElementHandler get_non_vol_element_handler(Coord x, Coord y, Coord z, CustomId id) const;
	
	inline static unsigned int is_empty_sphere(T &sphere_data) { return sphere_data == 0; }
};


// Devuelve el volumen parcial.
// Returns the accumulated volume
template <class T>
inline Volume Grillado<T>::obtener_vol_parcial() const { 
	return Volume(esferas)*vol_esfera-Volume(intersecciones)*vol_int; 
}

// Se calcula el volumen de las esferas.
// Calculates the volume of a sphere.
template <class T>
inline Volume Grillado<T>::calc_vol(Length radio) {
	return ((4.0f)*PI*radio*radio*radio/(3.0f)); }



// Se calcula el volumen de las intersecciones entre esferas.
// Calculates the volume of the intersection between neighbour spheres.
template <class T>
inline Volume Grillado<T>::calc_int(Length radio, Length dist) {
	return (2.0f)*PI*( ((2.0f)*radio*radio*radio/(3.0f)) - (dist*radio*radio/(2.0f)) + (dist*dist*dist/(24.0f)));
}


template <class T>
inline void Grillado<T>::reset_specific_sphere(T & sphere_data) const {
  sphere_data = 0;
}





//~ **********************************************************************************88








/*
    Spherical Grid Volume Approximation: An algorithm for approximating the 
    volume of a compound object by keeping track of its members' positions in space.
    Copyright (C) 2009, Rodrigo Castano, FuDePAN

    This file is part of Spherical Grid Volume Approximation.

    Spherical Grid Volume Approximation is free software: you can redistribute 
    it and/or modify it under the terms of the GNU General Public License as published 
    by the Free Software Foundation, either version 3 of the License, or
    (at your option) any later version.

    Spherical Grid Volume Approximation is distributed in the hope that it will be useful,
    but WITHOUT ANY WARRANTY; without even the implied warranty of
    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
    GNU General Public License for more details.

    You should have received a copy of the GNU General Public License
    along with this program.  If not, see <http://www.gnu.org/licenses/>.
*/


#include "grillado.h"

template <class T>
inline unsigned int Grillado<T>::show_quantity(T &sphere_data) const {
  return sphere_data;
}



// Se asigna memoria a la matriz y se inicializa con 0s.
// Para que no hayan intersecciones diagonales
// la suma de dos radios tiene que ser menor a la
// diagonal. Por Pitagoras, si x es la diagonal,
// x^2 = D^2 + D^2, luego buscamos que 2*R < sqrt(2)*D.
// Que es lo mismo que sqrt(2)*R < D.
template <class T>
Grillado<T>::Grillado(size_t M, size_t N, size_t Z, Length R, Length D) throw (radiusOrDistanceParamException, sizeParamException, bad_alloc)
	: v(M), w(N), z(Z), r(R), d(D), vol_esfera(calc_vol(R)), vol_int(calc_int(R, D)), intersecciones(0), esferas(0)
{
	const Length paramDiagonal = sqrt(2.0f)*R;
	const Length cota = (2.0f)*R;
	if(!((paramDiagonal)<D&&(D<cota))) {
		throw radiusOrDistanceParamException();
	}
	if(!(v>2||w>2||z>2)) {
		throw sizeParamException();
	}
	matriz = new T** [M];
	for(size_t m =0;m<M;m++) {
		matriz[m] = new T* [N];
		for(size_t n = 0;n<N;n++) {
			matriz[m][n] = new T[Z];
			for(size_t x=0;x<Z;x++) {
				reset_specific_sphere(matriz[m][n][x]);
			}
		}
	}
}


/* 	Se muestra el volumen de las esferas, el volumen de las intersecciones
 	y la cantidad de elementos en cada posicion del grillado.
*/
/* 	Shows the spheres' volume, the intersections' volume and the number of elements in each 
	position of the grid.
*/
template <class T>
void Grillado<T>::info_grillado() const {
	for(size_t m = 0;m<v;m++) {
		for(size_t n = 0;n<w;n++) {
			for(size_t Z = 0;Z<z;Z++) {
				cerr << "[" << show_quantity(matriz[m][n][Z]) << "]";
			}
			cerr << endl;
		}
		cerr << endl << endl;
	}
	cerr << "Sphere volume: " << vol_esfera << endl;
	cerr << "Intersection volume: " << vol_int << endl;
}
     
     
     
// Se libera la memoria correspondiente a la matriz
// All memory resources are released.
template <class T>
Grillado<T>::~Grillado<T>() { 
	for(size_t m = 0;m<v;m++) {
        	for(size_t n = 0;n<w;n++) {
                	delete [] matriz[m][n];
		}
		delete [] matriz[m];
	}
	delete [] matriz;
}






// Se quita un elemento de la esfera correspondiente a las coordenadas ingresadas
// como argumento.
// Removes an element from the sphere determined by coordinates x, y, z.





// Reduce el volumen parcial del grillado, se encapsula parte del comportamiento
// de sacar_esfera.
// Reduces the accumulated volume of the grid. It encapsulates part of the behaviour in
// sacar_esfera (removes a sphere)
template <class T>
void Grillado<T>::reducir_vol_parcial(GridCoord coord_x, GridCoord coord_y, GridCoord coord_z) {
	--esferas;
	intersecciones -= calcular_intersecciones(coord_x,coord_y,coord_z);
}


// Aumenta el volumen parcial del grillado, se encapsula parte del comportamiento
// de add_volumetric_element.
// Increases the accumulated volume of the grid. It encapsulates part of the behaviour in
// sacar_esfera (removes a sphere)
template <class T>
void Grillado<T>::aumentar_vol_parcial(GridCoord coord_x, GridCoord coord_y, GridCoord coord_z) {
	++esferas;
	intersecciones += calcular_intersecciones(coord_x,coord_y,coord_z);
}


// Asigna a coord_x, coord_y, coord_z los indices del grillado correspondientes a las 
// coordenadas x, y, z.
// Sets coord_x, coord_y, coord_z with the grid coordinates corresponding to x, y, z.
template <class T>
void Grillado<T>::calcular_coord(Coord x, Coord y, Coord Z, GridCoord &coord_x, GridCoord &coord_y, GridCoord &coord_z) const {
	coord_x = modulo(static_cast<int>(round(x / d )), static_cast<int>(v));
	coord_y = modulo(static_cast<int>(round(y / d )), static_cast<int>(w));
	coord_z = modulo(static_cast<int>(round(Z / d )), static_cast<int>(z));
}

/* 	Se calcula la cantidad de intersecciones de la esfera ubicada en el grillado en 
	coord_x, coord_y, coord_z con esferas vecinas que contengan algun elemento.
	Solo se tienen en cuenta las esferas cuyo centro pertenece a una recta paralela a 
	alguno de los ejes que ademas contenga al centro de la esfera cuyas intersecciones 
	se estan calculando.
	Es decir no se toman en cuenta intersecciones "diagonales". 
	Nota: estas intersecciones diagonales no pueden existir en la actual implementacion ya
	que las combinaciones de dist y radio que las generarian son consideradas
	invalidas por el constructor de la clase.
*/
/*
	Returns the number of intersections with direct neighbours of the sphere determined by coord_x, 
	coord_y, coord_z. 
*/
template <class T>
unsigned int Grillado<T>::calcular_intersecciones(GridCoord coord_x, GridCoord coord_y, GridCoord coord_z) const {
	//~ const GridCoord x_plus = modulo((static_cast<int>(coord_x) +1) , v);
	//~ const GridCoord x_minus = modulo((static_cast<int>(coord_x) - 1) , v);
	//~ const GridCoord y_minus = modulo((static_cast<int>(coord_y) - 1) , w);
	//~ const GridCoord y_plus = modulo((static_cast<int>(coord_y) + 1) , w);
	//~ const GridCoord z_plus = modulo((static_cast<int>(coord_z) + 1) , z);
	//~ const GridCoord z_minus = modulo((static_cast<int>(coord_z) - 1) , z);
	
	
	//~ unsigned int inters = !is_empty_sphere(matriz[coord_x][coord_y][z_plus]);
	//~ inters += !is_empty_sphere(matriz[coord_x][coord_y][z_minus]);
	//~ inters += !is_empty_sphere(matriz[x_plus][coord_y][coord_z]);
	//~ inters += !is_empty_sphere(matriz[x_minus][coord_y][coord_z]);
	//~ inters += !is_empty_sphere(matriz[coord_x][y_plus][coord_z]);
	//~ inters += !is_empty_sphere(matriz[coord_x][y_minus][coord_z]);
	
	const unsigned int inters = UTIL<Cross<1, 100, 100, 100 >::RES>::algorithm_list(matriz, coord_x, coord_y, coord_z );
	return inters;
}

// Sets every grid position to zero. Sets the partial volume to zero.
// Asigna cero a cada posicion del grillado y al volumen parcial.
template <class T>
void Grillado<T>::reset() {
	for(size_t m =0;m<v;m++) {
		for(size_t n = 0;n<w;n++) {
			for(size_t Z=0;Z<z;Z++) {
			  reset_specific_sphere(matriz[m][n][Z]);
			}
		}
	}
	esferas = 0;
	intersecciones = 0;
}



#endif // GRILLADO_H
